Effects of compression on the vibrational modes of marginally jammed solids 
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Glasses have a large excess of low-frequency vibrational modes in comparison with most crys- 
talline solids. We show that such a feature is a necessary consequence of the weak connectivity of 
the solid, and that the frequency of modes in excess is very sensitive to the pressure. We analyze 
in particular two systems whose density D{u!) of vibrational modes of angular frequency to display 
scaling behaviors with the packing fraction: (i) simulations of jammed packings of particles inter- 
acting through finite-range, purely repulsive potentials, comprised of weakly compressed spheres at 
zero temperature and (ii) a system with the same network of contacts, but where the force between 
any particles in contact (and therefore the total pressure) is set to zero. We account in the two cases 
for the observed a) convergence of -D(a;) toward a non-zero constant as — > 0, b) appearance of 
a low-frequency cutoff uj* , and c) power-law increase of ui* with compression. Differences between 
these two systems occur at lower frequency. The density of states of the modified system displays 
an abrupt plateau that appears at lj* , below which we expect the system to behave as a normal, 
continuous, elastic body. In the unmodified system, the pressure lowers the frequency of the modes 
in excess. The requirement of stability despite the destabilizing effect of pressure yields a lower 
bound on the number of extra contact per particle Sz: Sz > , which generalizes the Maxwell 
criterion for rigidity when pressure is present. This scaling behavior is observed in the simulations. 
We finally discuss how the cooling procedure can affect the microscopic structure and the density 
of normal modes. 

PACS numbers: 81.05.Rm, 82.70.-y, 83.80.Fg 
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I INTRODUCTION 



The responses to thermal or mechanical excitations 
and the transport properties of a solid depend on its low- 
frequency vibrations. In a continuous elastic body, trans- 
lational invariancc requires that the vibrational modes, 
or normal modes, are acoustic modes. This leads, no- 
tably, to the Debye behavior of the density of vibrational 
modes D{u!): in three dimensions, D{uj) ~ 0]. In 
amorphous solids structural disorder is present and this 
description breaks down at intermediate frequencies and 
length scales. In glasses, there is a large excess of low- 
frequency vibrations observed in neutron scattering ex- 
periments. This excess corresponds to a broad maximum 
in D{aj)/uj'^, called the "boson peak", which appears at 
frequencies on the order of a terahertz 0, that is typi- 
cally between one hundredth and one tenth of the Debye 
frequency. Above this frequency transport properties |5j 
arc strongly affected, as in silica where the thermal con- 
ductivity is several orders of magnitude smaller than in a 
crystal of identical composition [6j. Furthermore, amor- 
phous solids also exhibit force inhomogeneities that are 
common to granular materials and glasses 0. Thus, 
the resulting force propagation due to a localised pertur- 
bation differs considerably from the predictions of con- 
tinuum elasticity theory at short length scales. Recent 
simulations reported that this difference appears below a 



characteristic length scale which is directly related to the 
frequency of the boson peak 0|. Thus the excess- modes 
play an important role in several anomalous properties 
of amorphous solids. It was recently proposed that they 
also govern the physics at the liquid-glass transition |^. 

Despite their ubiquity, there is no accepted explana- 
tion of the underlying cause of these excess vibrational 
modes in glasses. A dramatic illustration of this excess 
was found in recent computer simulations of soft-spheres 
with repulsive, finite range potentials near the jamming 
transition . At this transition an amorphous solid loses 
both its bulk and shear moduli and becomes a liquid 0| . 
In a recent letter we showed how to calculate the 
density of states for weakly-connected amorphous solids, 
such as those near jamming, and showed that an excess 
density of vibrational states is a necessary feature of such 
systems. In this paper we use the method of [lH to pre- 
dict further consequences of compression. In particular, 
we explain why the coordination z must increase in a 
non-analytic way with applied pressure, as observed in 
the simulations US . 

Because there are no attractive forces and the temper- 
ature is zero in the soft-sphere simulations of O'Hern ct 
al- S S the pressure p = at the packing frac- 
tion at the jamming transition, cjjc- In 3-dimensions for 
mono-disperse spheres, (j)c ~ 0.64. The average number 
of contacting neighbors per particle, z, and the elastic 
moduli scale as functions of p. These simulations also 



reveal unexpected features in the density of vibrational 
modes, D{lo): (a) As shown in Fig. when the sys- 
tem is at the limit of marginal stability, as p ^ 0, D{uj) 
has a plateau extending down to zero frequency with no 
sign of the standard density of states normally ex- 
pected for a three-dimensional solid, (b) The plateau is 
progressively eroded at frequencies below a characteris- 
tic frequency w*, that increases with the pressure p (see 
Fig. 0}. (c) The value of D(ijj) of the plateau is un- 
affected by this compression, (d) At frequencies much 
lower than w*, D(lo) still increases much faster with uj 
than the quadratic Debye dependence. Similar behavior 
was also recently seen in models of tetrahedral covalent 
glasses earlier simulation of a Lennard- Jones 

glass had indicated that D{uj) increases at low frequen- 
cies when the coordination number, z, is lowered [l3j. In 
the following we aim to relate the density of vibrational 
modes of weakly-connected amorphous solids such as an 
assembly of finite-range, repulsive particles to their mi- 
croscopic structure. We start by briefly reviewing the 
results of 0| where the scaling properties of the density 
of states were derived. Then we use these results to com- 
pute the effect of applied stress. We show that it leads to 
a non-trivial constraint on the geometry of the contact 
network near the jamming threshold. The scaling behav- 
ior of the soft spheres near jamming furnish a stringent 
test to our predictions. 

At the center of our argument lies the concept of soft 
modes, or floppy modes. These are collective modes that 
conserve the distance, at first order, between any par- 
ticles in contact. They have been discussed in relation 
to various weakly-connected networks such as covalent 
glasses 0, 0, Alexander's models of soft solids [lij . 
models of static forces in granular packs [tR IT1 | and 
rigidity percolation models, see e.g. j20j. As we shall dis- 
cuss below, they are present when a system has too low 
a coordination number. As a consequence, as Maxwell 
showed |2]| . a system with a low average coordination 
number z has some soft modes and is therefore not rigid. 
There is a threshold value Zc where a system can become 
stable, such a state is called isostatic. As we shall dis- 
cuss this is the case at the jamming transition, if rattlers 
(particles with no contacts) are excluded. At this point 
there are no zero-frequency modes except for the trivial 
translation modes of the system as a whole. However, 
if any contact were to be removed, the frequency of one 
mode would go to zero , th at is, one soft mode would ap- 
pear. As we argued in this idea can be used to show 
that isostatic states have a constant density of states in 
any dimensions. When z > Zc, the system still behaves 
as an isostatic medium at a short enough length scale, 
which leads to the persistence of a plateau in the density 
of states at intermediate frequency. 

The second concept we use is at the heart of Alexan- 
der's discussion of soft solids 0|. In continuum elastic- 
ity the expansion of the energy for small displacements 
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FIG. 1: The density of vibrational states, D{ijj), vs angular 
frequency, uj, for the simulation of Ref. 0. 1024 spheres inter- 
acting with repulsive harmonic potentials were compressed in 
a periodic cubic box to packing fraction (p, slightly above the 
jamming threshold (l>c- Then the energy for arbitrary small 
displacements was calculated and the dynamical matrix in- 
ferred. The curve labeled a is at a relative packing fraction 
(j) — <f)c = 0.1. Proceeding to the left the curves have relative 
volume fractions 10^^, 10""^, 10~*, 10~*, respectively. 

contains a term proportional to the applied stress. It 
is responsible for the vibrations of strings and drum- 
heads and also for inelastic instabilities such as the buck- 
ling of thin rods. Alexander pointed out that this term 
also has a strong effect at a microscopic level in weakly- 
connected solids. For example, it confers rigidity to gels, 
even though these do not satisfy the Maxwell criterion 
for rigidity. We will show that while this term does not 
greatly affect the acoustic modes, it nevertheless strongly 
affects the soft modes. In a repulsive system of spheri- 
cal particles it lowers their frequency. We argue that 
this can dramatically change the density of states at low 
frequency, as confirmed by a comparison of simulations 
where the force in any contact is present, or set to zero. 
We show that these considerations also lead to a rela- 
tion between the excess connectivity 6z = z — Zc and the 
pressure, p. 

The paper is organized as follows. In the next section 
we write the expansion of the elastic energy that we use to 
derive the soft-modes equation. We then define and dis- 
cuss the isostatic case, and the nature of the soft modes 
that appear when contacts are removed in such a system. 
In the third part we compute the density of states when 
the effect of the applied stress on the vibrations is ne- 
glected. This approximation corresponds to a real physi- 
cal system: a network of relaxed springs. We use the soft 
modes in a variational argument to show that an isostatic 
state has a constant density of vibrational modes. We ex- 
tend this argument to include the case where where the 
coordination number increases above Zc, as is the case of 
the soft-sphere system under isotropic compression. We 
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show that such a system behaves as an isostatic one for 
length scales smaller than I* ~ Sz~^. This leads to a 
plateau in the density of states for frequency higher than 
LU* ^ Sz. At lower frequencies we expect the system to 
behave as a continuous medium with a Debye regime, 
which is consistent with our simulations. In the fourth 
part we study the effect of an applied pressure on D{uj). 
We show that although it does not affect the acoustic 
modes, it lowers the frequency of the soft modes. We 
give a simple scaling argument to evaluate this effect, 
and discuss its implication for the density of states. Inci- 
dentally this also furnishes an inequality between 5z and 
the pressure, which is also verified by the simulations, 
and which generalizes the Maxwell criterion for rigidity. 
Finally we discuss the influence of the cooling rate and 
the temperature history on the spatial structure and the 
density of states of the system, and conclude our work. 

II SOFT MODES AND ISOSTATICITY 
Energy expansion 

Following we consider N soft spheres packed into a 
spatially periodic cubic cell of side L at volume fraction 
(f). To describe the vibrations of a system of particles, 
wc expand the energy around equilibrium. For a central 
interaction V{r) one obtains: 

5E = Y, V'(r'^)dn, + \v'\rl])drl + 0{dr%) (1) 

ij 

where the sum is over all pairs of particles, r^' is the 
equilibrium distance between particles i and j. In order 
to get an expansion in terms of the particle displacements 
from equilibrium SRi...5Rn we use: 



dn, = [5R, - SRi).n,, + ^' ^ eg + 0{SR^) (2) 

where is the unit vector along the direction ij, and 
{6Rj —6Ri)^ indicates the projection of SRj — SRi on the 
plane orthogonal to Hij. When Eq. is used in Eq. J^l, 
the linear term in the displacement field disappears (the 
system is at equilibrium) and we obtain: 

+ lv"ir^J)[iSR, SR,).n,,r + 0(611') (4) 

In what follows we consider repulsive, finite-range "soft 
spheres" . For inter-particle distance r < cr, the parti- 



cles have non-zero mutual energy and are said to be in 
contact. They interact with the following potential: 

Vir) = -(l--Y (5) 

where a is the particle diameter and e a characteristic 
energy. For r > a the potential vanishes and particles do 
not interact. Henceforth we express all distances in units 
of (T, all energies in units of e, and all masses in units 
of the particle mass, to. In the following, we consider 
the harmonic case a = 1. In the discussion section we 
argue that these results can be extended for example to 
the case of Hertzian contacts |3| where a. = 5/2. In the 
harmonic case we have: 

^\ Y}^^^i - + 0(<Si?3) (6) 

where the sum is over all Nc contacts (ij). It is con- 
venient to express Eq. lO in matrix form, by defining 
the set of displacements 5Ri...SRn as a dA'^-componcnt 
vector |(5R). Then Eq. © can be written in the form 
SE = ((5R|A^|(5R). The corresponding matrix Ai is 
known as the dynamical matrix [3j. The dN eigenvec- 
tors of the dynamical matrix are the normal modes of 
the particle system, and its eigenvalues are the squared 
angular frequencies of these modes. 

The first term in Eq. © is proportional to the contact 
forces. In the rest of this paper we shall refer to this term 
as the stress term or transverse term. Near the jamming 
transition r^J 1 so that this term becomes arbitrarily 
small. We start by neglecting it, and we come back to 
its effects in the last section. This approximation corre- 
sponds to a real physical system where the soft spheres 
are replaced by point particles interacting with relaxed 
springs. We now have: 

M can be written as an by matrix whose elements 
are themselves tensors of rank d, the spatial dimension: 

Mij = -5(^j)nlj ® nlj + S.^j ^ n7; (E) nli (8) 

</> 

where (S^y) = 1 when i and j are in contact and the sum 
is taken on all the contacts I with i. 

Soft Modes 

If the system has too few contacts, A4 has a set of 
modes of vanishing restoring force and thus vanishing vi- 
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brational frequency. These are the soft modes mentioned 
in the introduction. For these soft modes the energy dE 
of Eq. Q must vanish; therefore they must satisfy the 
Nc constraint equations: 

{6Ri — 6Rj).nij = for all Nc contacts {ij) (9) 

This linear equation defines the vector space of displace- 
ment fields that conserve the distances at first order be- 
tween particles in contact. The particles can yield with- 
out restoring force if their displacements lie in this vector 
space. Eq. JH)) is purely geometrical and docs not depend 
on the interaction potential. Each equation restricts the 
dA^-dimensional space of |JR) by one dimension. In gen- 
eral, these dimensions arc independent, so that the num- 
ber of independent soft modes is dN — Nc- Of these, 
d{d+l)/2 modes arc dictated by the translational and ro- 
tational invariance of the energy function SElSR]. Apart 
from these, there are dN — Nc — d{d + l)/2 independent 
internal soft modes. 



Isostaticity 

There are no internal soft modes in a rigid solid. This 
is true when a system of repulsive spheres jams, when the 
rattlers (it i.e., particles without contacts) are removed. 
Therefore jammed states must satisfy Nc > dN — d{d + 
l)/2, which is the Maxwell criterion for rigidity. In fact at 
the jamming transition, if rattlers |29l | are removed, one 
can show this inequality becomes an equality [T^ll^ll^ . 
as was verified in Q. Such a system is called isostatic. 
The coordination number z is then Zc = 2Nc/N — > 2d. 

An isostatic system is marginally stable: if q contacts 
are cut, a space of soft modes of dimension q appears. For 
our coming argument we need to discuss the extended 
character of these modes. In general when only one con- 
tact (ij) is cut in an isostatic system, the corresponding 
soft mode is not localized near (ij). This arises from the 
non-locality of the isostatic condition that gives rise to 
the soft modes' as was confirmed in the isostatic simu- 
lations of Ref. [13, which observed that the amplitudes 
of the soft modes were spread over a nonzero fraction of 
the particles. When many contacts are severed, the ex- 
tended character of the soft modes that appear depends 
on the geometry of the region being cut. If this region 
is compact many of the soft modes are localized. For 
example cutting all the contacts inside a sphere totally 
disconnects each particle within the sphere. Most of the 
soft modes are then the individual translations of these 
particles and are not extended throughout the system. 

In what follows we will be particularly interested in the 
case where the region of the cut is a hyperplane as illus- 
trated in Fig. (0). In this situation occasionally particles 
in the vicinity of the hyperplane can be left with less than 
d contacts, so that trivial localized soft modes can also 
appear. However we expect that there is a non- vanishing 




FIG. 2: Illustration of the boundary contact removal process 
described in the text. Eighteen particles are confined in a 
square box of side L periodically continued horizontally and 
vertically. An isostatic packing requires 33 contacts in this 
two-dimensional system. An arbitrarily drawn vertical line 
divides the system. A contact is removed wherever the line 
separates the contact from the center of a particle. Twenty- 
eight small triangles mark the intact contacts; removed con- 
tacts are shown by the five white circles. 




FIG. 3: One soft mode in two dimensions for ~ 1000 
particles. Each particle is represented by a dot. The relative 
displacement of the soft mode is represented by a line seg- 
ment extending from the dot. The mode was created from a 
previously prepared isostatic configuration, periodic in both 
directions, following Q. 20 contacts along the vertical edges 
were then removed and the soft modes determined. The mode 
pictured is an arbitrary linear combination of these modes. 
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fraction q' of the total soft modes that are not locaHzed 
near the hyperplane, but rather extend over the entire 
system, hke the mode shown in Fig. ||2Jl. We shall define 
extended modes more precisely in the next section. 

Ill D(uj) IN A SYSTEM OF RELAXED SPRINGS 
Isostatic case 

Variational procedure 

We aim to show first that the density of states of an 
isostatic system does not vanish at zero frequency. D[lo) 
is the total number of modes per unit volume per unit 
frequency range. Therefore we have to show that there 
are at least on the order of ujL'^ normal modes with fre- 
quencies smaller than uj for any small w in a system of 
linear size L. As we justify later, if proven in a system of 
size L for lo ~ ui^ ~ this property can be extended 
to a larger range of lo indcpcndant of L. Therefore it is 
sufficient to show that they are of the order of L''^^ nor- 
mal modes with frequency of the order of instead of 
the order of one such mode in a continuous solid. 

To do so we use a variational argument: is a positive 
symmetric matrix. Therefore if a normalized mode has 
an energy 5E, we know that the lowest eigenmode has a 
frequency < V 5E. Such argument can be extended 
to a set of modes jSOj : if there are m orthonormal trial 
modes with energy 5E < oj^, then there are at least to/2 
eigcnmodcs with frequency smaller than \f2ujt. Therefore 
we are led to find of the order of L'^~^ trial orthonormal 
modes with energy of order . 

Trial modes 

Our procedure for identifying the lowest frequency 
modes resembles that used for an ordinary solid. An 
isolated block of solid has three soft modes that are sim- 
ply translations along the three co-ordinate axes. If the 
block is enclosed in a rigid container, translation is no 
longer a soft mode. However, one may recover the lowest- 
frequency, fundamental modes by making a smooth, si- 
nusoidal distortion of the original soft modes. We follow 
an analogous procedure to find the fundamental modes 
of our isostatic system. First we identify the soft modes 
associated with the boundary constraints by removing 
these constraints. Next we find a smooth, sinusoidal dis- 
tortion of these modes that allows us to restore these 
constraints. 

For concreteness we consider the three-dimensional cu- 
bical A^-particle system S of Ref . with periodic bound- 
ary conditions at the jamming threshold. We label the 
axes of the cube by x, y, z. S is isostatic, so that the re- 
moval of n contacts allows exactly n displacement modes 



with no restoring force. Consider for example the system 
S' built from S by removing the q ^ contacts cross- 
ing an arbitrary plane orthogonal to (ox) ; by convention 
at X = 0, see Fig. (01. 5', which has a free boundary 
condition instead of periodic ones a long (ox), contains 
a space of soft modes of dimension g[31|, instead of one 
such mode — the translation of the whole system — in a 
normal solid. As stated above, we suppose that a sub- 
space of dimension q' ^ I? of these soft modes contains 
only extended modes. We define the extension of a mode 
relative to the cut hyperplane in terms of the amplitudes 
of the mode at distance x from this hyperplane. Specifi- 
cally the extension e of a normalized mode |(5R) is defined 
by X;jSin^(^)(«l'5R)^ = e, where the notation (i|(5R) 
indicates the displacement of the particle i of the mode 
considered. For example, a uniform mode with (i|(5R) 
constant for all sites has e = 1/2 independent of L. On 
the other hand, if (i|(5R) = except for a site i adjacent 
to the cut hyperplane, the XijL ~ and e ~ We 
define the subspace of extended modes by setting a fixed 
threshold of extension eo of order 1 and thus including 
only soft modes /3 for which e/3 > eg. As we discussed 
in the last section, we expect that a fixed fraction of the 
soft modes remain extended as the system becomes large. 
Thus if g' is the dimension of the extended modes vec- 
tor space, we shall suppose that g'/g remains finite as 
L — > oo. The appendix presents our numerical evidence 
for this behavior. 

We now use the vector space of dimension g' ^ I? of 
extended soft modes of S' to build g' orthonormal trial 
modes of S of frequency of the order 1/L. Let us define 
I^R^) to be a normalized basis of this space, 1 < /3 < 
g'. These modes are not soft in the jammed system S 
since they deform the previous q contacts located near 
X = 0. Nevertheless a set of trial modes, |(5Rp, can 
still be formed by altering the soft modes so that they 
do not have an appreciable amplitude at the boundary 
where the contacts were severed. We seek to alter the soft 
mode to minimize the distortion at the severed contacts 
while minimizing the distortion elsewhere. Accordingly, 
for each soft mode /3 we define the corresponding trial- 
mode displacement (i|i5R*) to be: 

(z|<5R^) =C/3sin(^^)(*|JR^) (10) 

where the constants are introduced to normalize the 
modes. C/j depends of the spatial distribution of the 
mode /?. If for example, a highly localized mode has 
(i|(5R) = except for a site % adjacent to the cut plane, 
Cfi grows without bound as L ^ oo. In the case of 
extended modes CJ^"^ = Y.(i]) sin^(^)(jl'5R/3)^ = > 

— 1/2 

ep, and therefore Cp is bounded above by Cg . The 
sine factor suppresses the problematic gaps and overlaps 
at the q contacts near x ~ and x = L. The unit 
basis |(5R/3) can always be chosen such that the |(5Rp are 
orthogonal, simply because the modulation by a sine that 



relates the two sets is an invertible linear mapping in the 
subspace of extended modes. Furthermore one readily 
verifies that the energy of each |(5R^) is small, and that 
the sine modulation generates an energy of order 
as expected. Indeed we have from Eq. JT)): 

SE = C|^[(sin(^^)(*|5R^) -sin(^)0-|<5R^)).n,,f 

(11) 

Using Eq. (jSJ, and expanding the sine, one obtains: 
SE « C| ^cos2(^)^(n,, . e^nUmp) ■ n.,f (12) 

{ij} 

<e„-i(7r/L)2^(j|5R^)2(13) 

where e*r is the unit vector along (ox) , and where we used 
I cos I < 1. The sum on the contacts can be written as a 
sum on all the particles since only one index is present 
in each term. Using the normalization of the mode P 
and the fact that the coordination number of a sphere is 
bounded by a constant Zmax {zmax = 12 for 3 dimensional 
spheres [S^. ), one obtains: 



SE <e^\7r/Lfz^ax = ujl (14) 

We have found on the order of trial orthonormal 
modes of frequency bounded by w^, ^ ^/L, and we can 
apply the variational argument mentioned above: the av- 
erage density of states is bounded below by a constant 
below frequencies of the order lu^. In what follows, the 
trial modes introduced in Eq. (|10|l . which arc the soft 
modes modulated by a sine wave, shall be called "anoma- 
lous modes". 

To conclude, one may ask if this variational argument 
can be improved, for example by considering geometries 
of broken contacts different from the hyperplane surfaces 
we have considered so far. When contacts are cut to cre- 
ate a vector space of extended soft modes, the soft modes 
must be modulated with a function that vanishes where 
the contacts are broken in order to obtain trial modes 
of low energy. On the one hand, cutting many contacts 
increases the number of trial modes. On the other hand, 
if too many contacts are broken, the modulating func- 
tion must have many "nodes" where it vanishes. Con- 
sequently this function displays larger gradients and the 
energies of the trial modes increase. Cutting a surface 
(or many surfaces, as we shall discuss below) is the best 
compromise between these two opposit effects. Thus our 
argument gives a natural limit to the number of low- 
frequency states to be expected. 




5z 

FIG. 4: Scaling of uj* with the excess coordination number, 
Sz in the system with relaxed springs. The line has a slope 
one. 

Extension to a wider range of frequencies 

We may extend this argument to show that the bound 
on the average density of states extends to higher fre- 
quencies. If the cubic simulation box were now divided 
into sub-cubes of size L/m, each sub-cube must have 
a density of states equal to the same D{u)) as was de- 
rived above, but extending to frequencies on order of 
mujL. These subsystem modes must be present in the 
full system as well, therefore the bound on D{llj) extends 
to [0, ttilul]- We thus prove that the same bound on the 
average density of states holds down to sizes of the order 
of a few particles, corresponding to frequencies indepen- 
dent of L. We note that in d dimensions this argument 
may be repeated to yield a total number of modes, L'^~^, 
below a frequency lol ~ ^/L, thus yielding a limiting 
nonzero density of states in any dimension. 

We note that the trial modes of energy SE ^ l^^ that 
we introduce by cutting the full system into subsystems 
of size I are, by construction, localized to a distance scale 
I. Nevertheless we expect that these trial modes will hy- 
bridize with the trial modes of other, neighboring, sub- 
systems; the corresponding normal modes will therefore 
not to be localized to such short length scales. 



State with Sz > 

When the system is compressed and moves away from 
the jamming transition, the simulations show that the 
extra-coordination number Sz = z — Zc increases. In the 
simulation, the compression also increases the force at all 
contacts. However, in this section we will ignore these 
forces and focus our attention only on the contact net- 
work created by the compression. Any tension or com- 
pression in the contacts is removed. The effect on the 
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energy is to remove the first bracketed term in Eq. © 
above. We note that removing these forces, which must 
add to zero on each particle, does not disturb the equi- 
hbrium of tlie particles or create displacements. In this 
section we ignore the question of how bz depends on the 
degree of compression. We will return to this question in 
the next section. 

Compression causes AiVc — N5z/2 ~ L'^'Sz extra con- 
straints to appear in Eq. 0. Cutting the boundaries 
of the system, as we did above, relaxes q ~ L**"^ con- 
straints. For a large system, L'^Sz > L'^^^ and thus 
q < ANc- Thus Eq. (O is still over-constrained and 
there will be no soft modes in the system. However, 
as the systems become smaller, the excess number of 
const graints, ANc, diminishes; for L smaller than some 
I* ^ dz~^, q becomes larger than AA^cthe system is again 
under-constrained as was already noticed in |l7l |. This 
allows one to build low-frequency modes in subsystems 
smaller than I*. These modes appear above a cut-off 
frequency uj* ^ they are the "anomalous modes" 

that contribute to the flat plateau in D{lu) above lu*. In 
other words, anomalous modes with characteristic length 
smaller than I* are not affected very much by the extra 
contacts, and the density of states is unperturbed above 
a frequency uj* ^ Sz. This scaling is checked numerically 
in Fig.^ This prediction is in very good agreement with 
the data up to (5z w 2. 

At frequencies lower than uj* we expect the system to 
behave as a disordered, but not poorly-connected, elas- 
tic medium. We expect that the vibrational modes will 
be similar to the plane waves of a continuous elastic 
body. We refer to these modes as "acoustic modes". 
Thus we expect D{uj) at small uj to vary as uj'^^^c^'^, 
where c(dz) is the sound speed at the given compres- 
sion. This c may be inferred from the bulk and shear 
moduli measured in the simulations P, Q, > finds 
the transverse velocity ct ~ (Jz)^/^, and the longitudi- 
nal velocity c; ~ Sz'^ in both three- and two- dimensions. 
Thus at low frequency Diuj) is dominated by the trans- 
verse acoustic modes and a,t uj ~ uj* the acoustic density 
of states is uj'^~^c~'^ ~ 5z'^~^5z-'^/'^ - 5z'^/^-^. For a 
three-dimensional system the acoustic density of states 
should be dramatically smaller than the plateau density 
of states. Since there is no smooth connection between 
the two regimes we expect a sharp drop-off in D{uj) for 
UJ < UJ* . Such a drop-off is indeed observed, as seen for 
a three dimensional system in Fig. (O. In fact, because 
of the finite size of the simulation, no acoustic modes are 
apparent at w < cj* near the transition. 

The behavior of such systems near the jamming thresh- 
old thus depend on the frequency uj at which they are 
observed. For uj > uj* the system behaves as an isostatic 
system, and for uj < uj* it behaves as a continuous elastic 
medium. Since the transverse and the longitudinal ve- 
locities do not scale in the same way, the presence of a 
unique cross-over in frequency leads to the appearance of 




FIG. 5: Log-linear plot of the density of states for a 3-d system 
with N=1024 for tliree values oi cj> — 4>c in the soft-sphere 
system (dotted line) and the system where the applied stress 
term has been removed (solid line). 



two distinct length scales li and It, defined as li ^ ciuj*~^ 
and It ^ Ctuj*^^. These lengths correspond respectively 
to the wavelengths of the longitudinal and transverse 
acoustic modes at uj* . Note that since q ~ (5z°, one has 
li ^ I* . Interestingly, It ^ dz^^/"^ is the smallest system 
size for which acoustic modes exist. For smaller system 
sizes, the lowest frequency mode is not a plane wave, but 
is an anomalous mode. It can be observed numerically by 
considering the peak of the transverse structure factor at 
UJ* 0. 

Our argument ignores the spatial fluctuations of 5z. If 
these fluctuations were spatially uncorrelatcd they would 
be Gaussian upon coarse-graining: then the extra num- 
ber of contacts ANc in a subrcgion of size L would have 
fluctuations of order L''/^. The scaling of the contact 
number that appears in our description is ANc ~ L'^~^ 
and is therefore larger than these Gaussian fluctuations 
for d > 2. In other terms at the length scale I*, where 
soft modes appear, the fluctuations in the number of con- 
tacts inside the bulk arc negligible in comparison with the 
number of contacts at the surface. Therefore the anoma- 
lous modes are not sensitive to fluctuations in coordina- 
tion number in three dimensions near the transition. In 
the discussion section we will argue that there are spa- 
tial anti-correlations in z, so that fluctuations also do not 
affect the extended soft modes in two dimensions. 

Note that these arguments do not preclude the exis- 
tence of low-frequency localized modes that may appear 
in regions of small size I « I*, and that could be induced 
by very weak local coordination or specific arrangements 
of the particles. The presence of such modes would in- 
crease the density of states at low-frequency. There is no 
evidence for their presence in the simulations of 0. 
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IV EFFECTS OF THE APPLIED STRESS ON 
VIBRATIONAL MODES 

In this section we describe how the above simple de- 
scription of D{uj) is affected by the presence of an apphed 
stress. In general when a system of particles at equilib- 
rium is formed, there are forces between interacting parti- 
cles. For harmonic soft spheres it leads to a non- vanishing 
first term in Eq. ® iE(.y>(4' " ^W^j - SR,)^]\ 
where we used w 1. This term is (a) negative for 
repulsive particles (b) proportional to the transverse rel- 
ative displacement between particles in contact, (c) scales 
linearly with the pressure p and therefore vanishes at the 
jamming transition. The full dynamical matrix V can be 
written: 



V = M+M' (15) 

where A4' is written in tcnsorial notation in the footnote 
■ The spectrum of P has it a priori no simple relation 
with the spectrum of Ai. Because M' is much smaller 
than A4 near the transition, one can successfully use per- 
turbation theory for the bulk part of the normal modes 
of M ■ However perturbation theory fails at very low fre- 
quency, which is of most interest to us here. In this region 
the spectrum of A4 contains the acoustic modes and the 
"anomalous modes" forming the plateau. In what fol- 
lows we estimate the change of frequency induced by the 
applied stress on these modes. We show that the relative 
correction to the plane-wave frequencies is very small, 
whereas the frequency of the anomalous modes can be 
changed considerably. Finally we show that these con- 
siderations lead to a correction to the Maxwell rigidity 
criterion. 



Applied stress and acoustic modes 

Consider a plane wave of wave vector k. As the direc- 
tions fiij arc random, both the relative longitudinal and 
transverse displacements of this plane wave arc of the 
same order: [{SR^ - SRj)^]^ - [{SRi - 5Rj).nijf - . 
Consequently the relative correction AE/E induced by 
the applied stress term is very small: 
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FIG. 6: Sum of the transverse terms (full curve) e = 
1/2 [(5^i — SRi)^'^ and longitudinal terms (dotted 

curve) e = 1/2 {{^^j ~ ' ^sj]'^ for each mode of fre- 
quency Lo at the jamming threshold in 3-dimensions. The lon- 
gitudinal term is equal to the energy of the modes and vanish 
quadratically at frequency. The transverse term converges 
toward a constant different from 0. 

Applied stress and anomalous modes 

For anomalous modes the situation is very different: we 
expect the transverse relative displacements to be much 
larger than the longitudinal ones. Indeed soft modes were 
built by imposing zero longitudinal terms, but there were 
no constraints on the transverse ones. These are the de- 
grees of freedom that generate the large number of soft 
modes. The most simple assumption is that the rela- 
tive transverse displacements arc of the order of the dis- 
placements themselves, that is K'^-^j " ^Rd^V' ^ 
^ 2 

Ej 5Ri = 1 for the anomalous modes that appear above 
Lj* . This estimate can be checked numerically for an iso- 
static system where this sum is computed for all lo. The 
sum of the transverse relative displacements converges to 
a constant when w ^ as assumed, sec Fig. (0. 

Finally we can estimate the scaling of the correction in 
the energy A_B induced by the stress term on the anoma- 
lous modes: 

Ai? ~ - ^(1 - r^^WR, 6R,)^Y ^ -p (17) 

iv) 

which is an absolute correction, and which can be non- 
negligible in comparison with the energy E. 

Onset of the anomalous modes 



^ Ei:,)[m-6R,)-n,,]^ 

since r|^' — 1 is proportional to the pressure p, while the 
others factors remain constant as p 0, ~ p, and is 
thus arbitrarily small near the jamming threshold |34| . 



We can now estimate the lowest frequency of the 
anomalous modes. The modes that appear at lo* in 
the relaxed-spring system have an energy lowered by an 
amount on the order of —p in comparison to the origi- 
nal system. Applying the variational theorem of the last 
section to the collection of slow modes near w*, one finds 
that there must be slow normal modes with a lower en- 
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ergy. That is, the energy at which anomalous modes 
appear satisfies: 

iOAM^ < - A2P = AiSz^ - A2P (18) 

where Ai and A2 are two positive constants. Thus co- 
ordination and stress determine the onset of the excess- 
modes. 



Extended Maxwell criterion for the contact number 
under stress 

From this estimate we can readily obtain a relation 
between the coordination and pressure that guarantees 
the stabihty of a system. There should be no negative 
frequencies in a stable system, therefore loam > 0. Thus 
in an harmonic system the right hand side of Eq. H18() 
must be positive: 



Sz > Copl/" = 6Zrmn (19) 

where Co is a constant. This inequality, which must hold 
for any spatial dimension, indicates how a system must 
be connected to counterbalance the destabilizing effect of 
the pressure. A phase diagram of rigidity is represented 
in Fig. lO. When p = 0, the minimal coordination Zc 
corresponds to the isostatic state: this is the Maxwell 
criterion. As we said, for spherical particles Zc = 2d. 
When friction is present, one finds Zc = d + 1. When p > 
0, Eq. (O defimits the region of rigid systems: granular 
matter, emulsions lie above this line. As was shown by 
Alexander |l6l |. when p < 0, even systems with many 
fewer contacts than required by the Maxwell criterion are 
rigid. These systems contain many soft modes as defined 
in Eq. lO, but they are all stabilized by the positive 
bracketed term of Eq. ©. This is the case for example 
in a gel where polymers are stretched by the osmotic 
pressure of the solvent. Thus the network of reticulated 
polymers carries a negative pressure, which rigidifics the 
system and leads to a non-vanishing shear modulus. 

The relation H19() is verified in the simulations of Q 
where Sz ~ p^/^: the numerical results are in agreement 
with an equality of Eq. I|19|l . Furthermore it appears in 
Fig. (Ol that uj* » OJAM- In other words, the two op- 
posite effects of pressure on the vibrations, that is (i) the 
increase in the coordination and (ii) the addition of a 
negative term in the energy expansion, compensate |36|. 
In the discussion section we justify why the system of |3| 
is thus marginally stable even above the transition when 
p > 0, and we furnish examples of dynamics that lead to 
such features. 



Z 





. Granular Matter, 1 




emulsions^^.-'-^^ 








gels 

(Alexander) 





p=Q p 



FIG. 7: Phase diagram of rigidity in terms of coordination 
number and pressure. When p > 0, the line separating the 
stable and unstable regions is defined by Eq. When the 

pressure is negative, any connected system can be rigid. 

DISCUSSION 
Stronger constraints on Sz 

The simulations of show Sz ~ p^^^, thus poten- 
tially saturating the bound of Eq. (20), so that there 
are excess modes extending to frequencies much lower 
than Lo* . Here we furnish an example of dynamics that 
lead to such a situation. Consider an initial condition 
where forces are balanced on every particle, but such 
that the inequality (|19() is not satisfied. Consequently, 
this system is not stable: infinitesimal fluctuations make 
the system relax with the collapse of unstable modes. 
Such dynamics were described by Alexander in |l(il | as 
structural buckling events: these events are induced by 
a positive stress as occurs in the buckling of a rod; but 
here they take place in the bulk of an amorphous solid. 
These events a priori create both new contacts and de- 
crease the pressure. When the bound of lfTn|) is reached, 
there are no more unstable modes. If the temperature 
is zero, the dynamics stop. Consequently one obtains a 
system where Eq. ()19|l is an equality. Therefore (i) this 
system is weakly connected (ii) loam ~ 0, so that there 
arc anomalous modes very different from acoustic modes 
extending to zero frequency. 

In the simulations of Ref. 0] the relaxation proceeds in 
the following way. The system is initially in equilibrium 
at a high temperature. Then it is quenched instanta- 
neously to zero temperature. At short time scales the 
dynamics that follows is dominated by the relaxation of 
the stable, high frequency modes. The main effect is to 
restore approximately force balance on every particle. At 
this point, if the inequality H19|) is satisfied, then the dy- 
namics stop. However, if it is not satisfied then we are in 
the situation of buckling as described above. The pres- 
sure and co-ordination number continue to change until 
the last unstable mode has been stabilized. At this point 
the bound of Eq. (20) is marginally satisfied, and there 
is no driving force for further relaxation. 

The data of |2| were obtained by gradually decreasing 
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the pressure from this initial state of zero temperature 
and nonzero pressure. The reduction of pressure causes 
some contacts to open. The opening of these contacts 
tends to destabilize normal modes and reduce their fre- 
quencies, while the reduction in pressure tends to stabi- 
lize them. If the particles simply spread apart affinely, 
the destabilizing effect would be expected to dominate 
[Sfij l . Thus we argue that some incremental buckling must 
occur as in the initial temperature quench. The buckling 
increases the contact number and decreases the pressure 
until marginal stability is achieved, so that the inequal- 
ity of Eq. H19() is marginally satisfied as the pressure de- 
creases. 

It is interesting to discuss further which systems, and 
which procedures, can have marginal stable states. In 
repulsive short-range systems, we expect the situation of 
marginal stability that follows an infinite cooling rate to 
take place for a domain of the parameters of initial condi- 
tions ((/), T), located at high temperature and low density. 
This domain might stop at a finite even when the tem- 
perature is infinite. This is true for particles interacting 
with a gaussian potential, as was shown in simulations 
and theoretical analysis on Euclidian random matrices 
|23j . There most of the unstable modes vanish at a finite 
(j) despite T = oo. 

When the cooling rate is finite we expect that the re- 
laxation does not stop when all the modes are stable. For 
example activated events can in principle lead to the col- 
lapse of anomalous modes. These events a priori increase 
the connectivity and decrease the pressure further than 
the bound of Eq. (|19|l . leading to ujam > 0. Interest- 
ingly hyper-quenched mineral glasses show a much larger 
number of excess modes |24l | in comparison with normally 
cooled glasses, whereas annealed polymeric glasses show 
a diminished excess ofsuch modes l26l . 



Fluctuation of connectivity 

In Section III we argued that although fiuctuations 
were negligible in three dimensions near the transition, 
they might have effects in two dimensions. Here we ad- 
dress these possible effects. In other systems such as fer- 
romagnets the spatial fluctuations of eg. magnetic cou- 
pling lead to noticeable and sometimes striking effects 
[23 . It is natural to ask what analogous effects might be 
induced by fluctuations of contact number z in the sys- 
tem of in two dimensions. In a system with uncorre- 
lated contacts in two dimensions, the modes in the weakly 
connected regions would be softer than those in the more 
strongly connected regions. This would be expected to 
create modes with frequencies below the nominal lu* in 
the system of relaxed springs of Section III, and would 
imply the presence of modes below the bound loam in the 
original system. We argue that in the simulations of 
such fluctuations are negligible. It is important to note 



that z is not an uncorrelated random variable like the 
magnetic coupling mentioned above. For example, if z 
were an uncorrelated random variable, one would expect 
that the fluctuations in the total number of contacts Nz 
would be of order However, at the jamming thresh- 
old the system is isostatic and the number of contacts is 
precisely 2dN; there are no fluctuations. More generally, 
the bound on Sz in the last section applies not merely 
on average but rigorously to any subregion of size larger 
than since all sub- regions experience the same pres- 
sure p. Furthermore for a marginally stable system such 
as those of 0, the average 6z is given by the bound of 
Eq. (|19|) . These two facts imply that the fluctuations of 
contact number are small. Therefore the conclusions of 
the last section remain valid, and no anomalous modes 
below the bound uj* are expected in the relaxed-spring 
system due to fluctuations of z. 

Extension to non-harmonic contacts 

In the previous sections we considered harmonic inter- 
actions. Here we discuss the generalization of our argu- 
ment to other potentials. Ref. explored several other 
interactions, notably the Hertzian interaction potential 
describing the compressive energy of two elastic spheres. 
It corresponds to a = 5/2 in Eq. (jSJ. Ref. |2| observed a 
plateau in the density of states whose height Dq scales as 
They also observed a cutoff frequency uj* varying 
as p^/^. In the Herztian case the quadratic energy of Eq. 
(7) becomes: 

SE^^ E(l - ^v)'^'li^Rj - SR^).n^J? (20) 

The new factor (1 — rijY^'^ amounts to a spring con- 
stant kij that depends on compression. The contact force 
fij = dSE/drij evidently varies as (1 — rijf'^'^ ■ In what 
follows we neglect the fluctuations that exist between the 
contacts. This treatment is sufficient to recover the scal- 
ing results of . 

The new factor (1 — Hj)^^^ rescales the energy. To 
account for this overall effect, we replace (1 — rij)^!"^ by 
its average ((1 — ri^f-^"^)- Expressed in terms of contact 

1 /3 

forces, this factor is proportional to (/^^ ). Replacing 
by its average, the factor becomes {fijY^^ ■ This average 
is related to the pressure p, via p « {fij)- Thus, in this 
approximation the overall effect is to rescale the energy 
by a factor k{p) ~ p^/'^. 

SE=Mj2[(^SR,~SR.).ri,,r (21) 

Apart from this prefactor, the energy and the dynam- 
ical matrix are identical to the harmonic case treated 
above. Each normal mode frequency gains a factor 
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^ In the harmonic case the crossover fre- 

quency foUows uj* ^ Sz. In the Hertzian case, it gains 
the same factor k^^^, so that ui* ^ k^/^5z. The bound 
on the lowest-frequency anomalous modes uJj^m still has 
the form 

tOAM^ < io*^ - A2P (22) 

For a marginally stable system we still have ujam ~ 0, 
which leads to an unaltered relationship between lu* and 
p: UJ* ^ p^l"^ . Comparing with our previous estimate of 
Lo* we find 8z ^ p^^^ . Furthermore, the plateau density 
of states Dq has the dimension of an inverse frequency 
and thus gains a factor Since the harmonic Dq had 

no dependence on p, the Hertzian (p) also should vary 
as p~^^^. The scaling behaviors seen in j3| agree with 
these expectations. These arguments may be applied to 
general values of the interaction exponent a. 

Additional effects could in principle alter the low- 
frequency modes in the Hertzian case. When harmonic 
springs are replaced by Hertzian ones, different contacts 
have different stiffness. This effect should be quantified 
in order to gain a more detailed understanding of the 
Hertzian case |2^ . 

CONCLUSION 



The anomalous modes offer a new approach for un- 
derstanding response and transport in weakly-connected 
mechanical systems. Knowledge of the statistical proper- 
ties of the anomalous and soft modes defined by Eq. 
is necessary to predict the acoustic and thermal trans- 
port of e.g. the simulated system of Q]. These modes 
also provide an alternative view of relaxation, via buck- 
ling of compressed anomalous modes. This buckling pic- 
ture provides an intriguing contrast to the local cage- 
escape picture commonly used to describe these relax- 
ations. This approach may be more broadly applicable 
to glasses. Such glasses exhibit a large excess of low fre- 
quency modes like our marginally jammed system; this 
suggests that they behave like weakly-connected mechan- 
ical systems at short length scales. This resemblance 
raises the hope that our approach may explain anoma- 
lous phenomena in the transport, response, relaxation 
and aging of structural glasses. 
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edge the support of the CFR fellowship, MRSEC DMR- 
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In this paper we computed some vibrational properties 
of weakly-connected amorphous systems. We introduced 
a frequency scale cu* above which such systems do not 
behave as continuous elastic media, but as isostatic sys- 
tems that are marginally stable. At frequencies lower 
than UJ* , the system acts as a continuous solid, with low- 
frequency, acoustic, phonon modes. At frequencies above 
UJ* the vibrational modes correspond to the anomalous 
modes that constitute the vibrations of isostatic states. 
As we showed in Section „ these anomalous modes arc 
built from the soft modes that appear in subsystems with 
free-boundary conditions. The anomalous modes that 
appear at w* are characterized by a length scale /*. In- 
terestingly I* does not appear directly in any correlation 
function of the static structure, but rather in the response 
fucntions of such systems. I* can be much larger than the 
particle size and varies with the coordination number, as 
I* ^ uj*~^ ^ Sz~^. Secondly we computed the effect 
of applied stress on these anomalous modes. In a re- 
pulsive system the stress has a strong effect that lowers 
the frequency of the anomalous modes. Imposing that 
such modes are stable leads to a generalization of the 
Maxwell criterion: the coordination must increase non- 
analytically with pressure to compensate for the desta- 
bilizing effect of compression. Finally, we discussed the 
"structural buckling" that occurs when anomalous modes 
collapse. We use this concept to justify the marginal sta- 
bility that follows hyper- quenches in the simulations of 
repulsive, short-range systems of 



APPENDIX: SPATIAL DISTRIBUTION OF THE 
SOFT MODES 

In our argument we have assumed that when q ~ L"^^^ 
contacts were cut along a hyperplane in an isostatic sys- 
tem, there was a vector space of dimension q' = aq which 
contains only extended modes, and that a does not van- 
ish when L 00. A normalized mode |(5R) was said to 
be extended if ^ ■ sin^(^)(i|5R)^ > Cq, where eg is a 
constant, and does not depend on L. Here we show how 
to chose eg > so that there is a non-vanishing fraction 
of extended soft modes. We build the vector space of ex- 
tended soft modes and furnish a bound to its dimension. 

Let us consider the linear mapping Q which assigns to a 
displacement field |(5R) the displacement field (i|fJ(5R) = 
sin^ {xi'K/L)SRi. For any soft mode |(5R/3) one can 
consider the positive number ap = {SUplQlSHfs) = 
sin^(xi7r/L)(5_R^^. We build the vector space of ex- 
tended modes by recurrence: at each step we compute 
the a/3 for the normalized soft modes, and we eliminate 
the soft mode with the minimum a/3. We then repeat 
this procedure in the vector space orthogonal to the soft 
modes eliminated. We stop the procedure when a/3 > eg 
for all the soft modes /3 left. Then all the modes left 
are extended according to our definition. We just have 
to show that one can choose eo > such that when this 
procedure stops, there are q' modes left, with q' > aq, 
with a > 0. In order to show that, we introduce the 
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smaller than 1. Thus X]^=i '^P < {q ^ q')so + q' ■ Com- 
bining, we have {q — q')eQ + q' > qc/2, or 
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FIG. 8: The overlap function f{x) as defined in Eq. I|23^ for 
different system sizes in three dimensions. Soft modes were 
created from isostatic configurations as described in Fig. 



following overlap function: 



f{x)dx^q-' E ['^^^.^ 



(23) 



/3— !,...,(? Xi^[x,x+dx\ 



The sum is taken on an orthonormal basis of soft modes 
/3 and on all the particles whose position has a coordi- 
nate Xi e [cCjX -|- dx]. f{x) is the trace of a projection, 
and is therefore independent of the orthonormal basis 
considered. f{x) describes the spatial distribution of the 
amplitude of the soft modes. The |(5R/3) are normalised 
and therefore: 



f{x)dx 



(24) 



g' >(7(c/2 - eo)/(l - eo) 



(26) 



Evidently for all L, q' remains a nonvanishing fraction 
of q for any fixed threshold eo such that eo < c/2, as 
claimed. 
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[32] In a polydisperse system Zmax could a priori be larger. 
Nevertheless Eq. I|13|l is a sum on every contact where 
the displacement of only one of the two particles appears 
in each term of the sum. The corresponding particle can 
be chosen arbitrarily. It is convenient to choose the small- 
est particle of each contact. Thus when this sum on ev- 
ery contact is written as a sum on every particle to ob- 
tain Eq. (I14II . the constant Zmax still corresponds to the 
monodisperse particle cannot have more con- 

tacts than that with particles larger than itself. 

[33] in three dimension we have M'ij = — ^2r '^-' [^(iJ) ('^j ® 
rriij + kij ®kij) + 5i,j J2<:i^i''^ij ®rnij +kij ®kij), where 
{riij ,rnij ,kij) is an orthonormal basis. 

[34] In disordered systems the acoustic modes are not exact 
acoustic modes, see e.g. recent simulations in Lennard- 



Jones systems Q- Therefore the correction due to the 
applied pressure on the frequency of the acoustic modes 
might be larger than what is expected from Eq. 1161 . 
Nevertheless, we still expect the effect of applied pressure 
to be much smaller on the acoustic modes than on the 
anomalous modes. 

[36] Assuming an exact compensation of these two terms lead 
to ljam = in an infinite size system. In Fig. ujsm is 
slightly different from zero. We believe that it is simply 
a finite size effect. 

[36] When the particle's radius decreases by an amount e, 
a certain fraction e ~ 5(l)e of contacts opens, where 
g{r) denotes the radial distribution function. For har- 
monic particles we expect g{l) ^ {(f) — (/)c)~^, as we jus- 
tify at the end of the paragraph. Thus using the fact 
that p ^ {(j} — (jjc) for harmonic particles [J, one ob- 
tains e ^ . On the other hand, the affine decompression 
lowers the pressure by an amount 5p S(j> ^ e. Thus 
the fraction / of contacts that the system can afford to 
lose while staying stable follows, according to Eq. 1191 : 
/ = i^%^5p ~ p-i/^Jp ^ « e. Therefore / « e 
as claimed. Thus if an affine reduction of packing fraction 
e is imposed, far too many contacts open and the system 
is unstable. To conclude we justify g{l) ^ [(f) — 0c)~^. 
This is related to well-known empirical facts of the force 
distribution: one has P{F)dF ~ g{r)dr. For harmonic 
particles dr ~ dF and therefore P{F) ~ gix)- When 
rescaled by {F) ~ Pi^P{F) converges to a master curve 
with P{F = 0) / la,Il. This implies that g{l) ~ p'^ ~ 



